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ABSTRACT 


USinwmoncewmmeve ll barotropre ocean model, driven by sur- 
PimVvilas, 2eriiiec Gifterence form Of the vorticity equation 
Poomantegrated over 210 days of simulated time. The solutions 
using constant coefficients of lateral eddy viscosity were 
moloanred With those using variable coefficients derived from 
PMourOphy cascade and energy cascade. Using a constant eddy 
ios ity Coefficient of rather low magnitude produces a large 
Dpicnae Computational oscillation which fills the entire 
wipe An Order of magnitude larger coefficient produces a 
feaganially satisfactory solution, where the western boundary 
Mperene was rather welll represented, but a moderate computa- 
meonal oscillation was still evident. By increasing the co- 
efficient yet another order of magnitude, the computational 
Peeillation is negligible, but the solution in the ocean in- 
memnor 1S unrealistically damped. An accurate physical and 
iemerrecal depiction of both the ocean interior and western 
boundary with no computational oscillation was achieved by 


meege Cither of the two forms of non linear eddy viscosity. 
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I. INTRODUCTION 


HiiSmeaicdy smart eo a broader €ffort to deve lopathe 
Gapapility ot making large scale oceanographic predictions 
Sia dynamical basis. The large-scale and meso-scale thermal 
femucture of the ocean 1S very dynamze and is related to most 
iarine and atmospheric processes. Some of the most obvious 
relationships to sea surface temperature (SST) are ocean 
meones, Clrculation, currents, and Sea state; and atmospheric 
meniperaeures Clreubatrons, @and witid velocity. The effect that 
these environmental phenomena can have on naval and maritime 
Operations is well» known. 

The significance of the long range effect of SST anomalies 
Gumweather patterns has received increasing scientific aware- 
Mewes On 20 August 1974, J. Namias, at a NORPAX Conference 
mampeic U.S. Naval Postgraduate School, reported on the results 
Sean empirical ocean/atmosphere prediction model. With a 
mewbedge of SST anomalies in the spring of 1974, Namias 
memived fields of atmospheric temperature anomalies and cir- 
culation for the following summer. Namias predicted in May 
fae a comparatively severe drought would occur over the 
meawestemn United States during the summer of 1974. Later 
Peet sceyvertirced his forecast. 

One worethiecmancadcmor dittlculeyeam adequately representing 
MiceCl Vel HOnmalmdmanonalies in @ 2inite-crid ocean circula- 


tion model, and subsequently satisfying and integrating the 





Pecos wOreMOeron, 1S the representation, of internal 
PeletVronalwecmececs within the ocean fluid. 

cme lincartheory, Takamew{1974]) showed that if the 
EGeriteleme OLeecddy ViSCOSity 1S too small, a computational 
Tioveeerene iiespace results from not resolvina the western 
Dome, sclrrent.. this computational oscr1llation f£111s the 
pierre basin awd comtaminatess the Solution. Onethe other 
Min@eniiwtne "COC flGlent 1S increadsed, the open ocean soln- 
len 15 too viscous. To handle this problem, Takano showed 
midtmeene USe Of UpStream dinterencine im the beta term of 
mie VOrtureity equation alllows a smaller coefficient than per- 
ipeercduby linear theory. Hewever, this method unfortunately 
ireoauces excessive damping of time dependent motion. In 
madditr1on, this mathematical scheme 1s not representative of 
mee physical motion in the ocean. <A better approach, from 
Emon ysitcal point Of view, 1S the use of non linear eddy vis- 
Beatty. the friction force which arises using nom linear 
eer ViISCOsity 15 Extremely sensitive to the wcale of the 
Meelon, and therefore will be relatively small in the oceans' 
meerlon where the scales are comparatively large, and will 
bemrelatively large in the western boundary region where the 
scales are comparatively small. The comparatively large 
dissipation in the western boundary region will een ihe 
wrirenme road enough Lo be resolved by the grid and thereby 
Meech EP eCNewLOrMatlon Of any Computational oscillation. 

In) #968, Leith examined two dimensional turbulence 


advection and derived non linear coefficients of eddy 


10 





VaroeCs lt) trom tiemedscade Of enstrophy from large scale to 
Miieilececale mMotlommeeliineniS case, the coefficient of eddy 
MEscOSTLy 1Sapnoperulonal to the magnitude of the horizontal 
Sraalene Of VOrtlel cy. 

AUCthetas Mees cetonry mon linear procedure was introduced 
maonmacorinsky | 6905] an which an estimate of the energy cas - 
@iact rate 15 Mage from the fluid @eformation itself. In this 
ao Cec cmimlbetcniGr caddy VISCOSILY 15 proportional to 
Grewaesolute value of the deformation. 

iiemDUGhOsecuot stiISh thesis: is to introduce non linear 
MOcrteclents Of lateral eddy viscosity that are not only 
fmeasonably representative of actual physical conditions, but 
Bees ceSane COctfIGClents must retain their reasonable repre- 
meiearcaon when used in finite grid numerical ocean circulation 
Mmedels., In addition, with these coefficients the numerical 
Peo@ecle must remain computationally stable in long term inte- 
fegeron. Integrations using both of the above forms of non 
migear Cddy viscosity are €xamined and compared with numerical 


aa analytic results using constant coefficients. 


ied 





Piet MAT ew eA oT A DEMENT VOR THE PROBLEM 


Pee FORM OF THE VORTICITY EQUATEON 

ice contin, MeQuathoOlNerO. aulomovcneoous fluid and the 
WW inear Calaerons Of horizontal motion, cross differenti - 
a ecCOecriiiiatcet lcm UnhessineNrcermn, Lor the basis forthe 


Womelecilty Caulation: 


Er 


‘s 
=—yt) oe ~~) 


0 rt 
dt + vB + fv-V = 9oxl¥y + oo) - TyCFy * DoH 


IX 


Or 


Q? 
ct 


3 
*V-ve pe ee Paty a 5 Fx Fe 
(i= 1) 
where g = as Pee keioeconstantcullm (ll-1), {'x, ty) aS the 
Surface stress, H is the constant basin depth, and all other 
symbols have their usual Neate. — BOttomeust:ess was meolected- 

Wtcumerc tlOWmnROYTCCSmare TeMrescnted 

a cy 2) 

Fy = y- (Ay v) | (11-2) 
where A is the coefficient of lateral eddy viscosity and (u,v) 
/memene Vertically avyerageasvelocitvemas discussedebelow. 

Iiesmodel isedesiened with a “rigid surface,” in order 
Memtilter gravity waves from the system: 

w(o) = 0. 
ierts also Mecessary that there be no vertical motion on the 
fbat ocean floor: 

w(-H) = 0. 


Peealsemw els szZero et the top and bottom, the divergence of 


eZ 





Giemuertically averaced Gurrent, V, is zero. 


YooV = 0. 
Therefore V can be defined by a streamfunction, pt , 

fe oe! 

oy 

ee! du 

\X. 


Integrating (II-1) from top to bottom and using y+ V = 0, 
the final form of the vertically integrated vorticity equation 


was written 








oe Ds ee ers ay. 
ee Py Ge VY + Fy * coe 
Ula aor, = ex 
merece OPE on) nes) 


lpm ON Wl1-3)esalso includes ythe assumption that 
Wey =V«9v, etc. 
Hitcmwiecessecenon will diseuss the formulation of the 


merction terms, F, and F VirechiadencmamuMnOlmele Non. Linear 


y? 
Saad Viscosity coefficient, A. 
fee ENS TROPHY CASCADE FORM OF NON LINEAR EDDY VISCOSITY 
iiemimirstaietnod Of Peneratimea non linear coerficients 
Merecddy ViScCOsity 15 through the theory of two dimensional 
Mipmlence ite witGh EMStropny amd Kinetic energy are con- 
eemved py the advective terms. Using these principles, 
iMeadiecnman |l96/) and Leith [1968] derived the relation 
E(k) « n¥k73 C= ay 
where n is the assumed constant cascading rate of mean squared 
Mommbeiiny saldwh(m@jeas tne Kimetie energy in wave number k. 


The eddy viscosity which causes the dissipation is also assumed 


iS 





Peete rlon OF yoanmd k;, and by dimensional analysis 

A =a 7¥3k-? (11-5) 
Miicmeme LS a CONS tan t. 

Oiesestimave Grenjpmade locally as a £unction Of Sur- 
rounding data was made by Leith 

mn = NC We ° (CPE) SN el CaIS6) 
Piece GC is Eine local) Vvereical component of vorticity. 
pms theucime (li-6) into (11-5) and solving for A leads to 

A= (a ¥/k)?|Vel. 
Assuming that the wave number ky = 27/d, where d is the grid 
suze, lies in the inertial subrange, the final form of non 
iimicear €ddy viscosity for enstrophy cascade becomes 


= dy) 3 _ 
I (0 )/ xa) vel. | (II-7) 


omoting the Minimum viscosity coefficient by A,, the 


viscosity in the model was written, 


hea + yivcld- (11-8) 
A is a maximum when |Vc| is maximum. Defining Amax = mg; 


then (II-8) becomes 
Amax = MAg = Ao +¥1VSlmax 4° Giir=9)) 


Sormevying for y 


ema 1) : 
Y Sivas a? 
Barone Gielen 


ma 


ive Guat tyaiVyo| max Was estimated by 


sie moron 059 cn Dcec-! (II-11) 


Ticino red > 300ekm leads) to a value for. of 
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Ao Un-1) 


Y = 
Je ee 


Gane) 


PBs tmiiting thes wanto (1-8) "ives the final cquation for 


Micmenstnophy Cascade formmeo: mon linear eddy viscosity for 


Cars study: 
A = AQ(1 + Gap ively = Ag(1+ Ged) [Ve Gliees 
Vol max Ne Wy 


ieCii-13), ¥_ was taken from linear theory [Munk, 1950], 


ax 
and mM was considered an adjustable parameter. 
MitiieenisenOnelanear form of eddy viscosity, thie Wen tients, 


Mare Jeaerr 1 Gtlons force terms in the x end y direetions become 


trom (11-2) 


Be = eZ) yy 2 Ao) (ime, 
ox aX Jy dy 

eee ye bts 

y = 5x Aa? ih yy? ( 


em Chul Co weconmonent Orwene Curl of the friction force 


becomes 
Git, Se es | CU = 
ox oy 


eee x INETIC ENERGY CASCADE FORM OF NON LINEAR EDDY VISCOSITY 
From diminsional arguments, Leith [1968] showed that in 

aaoace smc stonal turbulence, 1£ the enérgy cascade rate from 

maee scale £6 Small Scale 1s proportional only to kinetic 


energy dissipation (€) and wave number (k), then 

Ck) 2 en eae ee (11-17) 
Mremeiimaddiiien to this equation, it 1S assumed that eddy 
Viscosity, which produces the dissipation, is also a function 


of ¢ and k only, then by dimensional arguments 
iS 





A =a, ek ¥ (11-18) 
where k lies in an inertial subrange. Assuming the dissipa- 
tion rate is constant, then k°¥3.d¥3. This quasi-linear 
Poh SCOSlLY ts udependent only on grad size, and hence 
Piltide 1: the grid size ts latitude dependent. But more 
realistically, dissipation is not constant nor independent 
of motion or grid size. Smagorinsky {1963} assumes a local 
male Of dissipation: 

2 = NDE (11-19) 
where |D| is the magnitude of the deformation tensor. In 


this case (II-18) becomes 


R= fp 12 10) Ge (11-20) 
Siteewe deformatron tensor |D| = / D7 +D,* , the shearing 
deformation is Dy = 9% +.94. and stretching deformation is 

| . dx ay 
De 


The Smagorinsky equation takes the form similar to 
wer-S) for this study 


he Net lid = (11-21) 


Reining AA WAG: equation (11-21) becomes 
me = Ne oP SPD Eva cle Li = 22) 
Solving for y 


y = Aotm-1) . (11-23) 


| DI max dé 


Usanie (1l-235) am Tl-21); the final energy cascade form of 


VicmioMmemimicanrucddyveviScosity becomes 
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- Poe = oD 
A= A, [1 inet tt (Gua) 


it mriOrmiaore (nest TiGthoimtorce Principally used im the 
model for the kinetic energy cascade case was of the form 
Tet ieandas(ll-15).9 Experiments following the frietion 


force Of omagorinsky [1963] were also used for comparison: 


Fx = (A D+) + -- (A De) (II-25) 
9x ay 
Een = 8 (A 1 = 02 AD IT- 26 


imo rn cases the curl of the friction force term was (11-16). 


i, 





Cw ieee DIP PERENCE EOUATTIONS 





Peed YorCAl CHARACTERISTICS OR@THE MODEL 

hie ued CiCcoOns ISted. Of a One - lever barotropie ocean a 42 
Square basin of length, L = 9600 km; of breadth, B = 9600 kn; 
imiced that botlom Of constant depth, H = 2 km. A portion of 
Mmremstaggened orld 1s shown in Figure 1. The total grid 
memmts (excluding the boundary buffer) are 33 x 33, and the 
distance between adjacent grid points is 

X = Y = d = 300 km. C= 1h) 
Bieevalue chosen for gp Corresponds to a grid centered at 
Byeo” N. 

icwmere Comednonnd CONSIStS Primarily of: 1) the inter- 
PeeeLOnsS sor Primcipal grid points (+ ), where are defined 
fijemstreaneunction (¥ ), vorticity (7), deformation (D), 
and the coefficient of eddy viscosity (Ap) generated from 
deformation; and 2) the grid centers (0), where are defined 
Pemocity (U,V), the gradient of vorticity (V.z), the friction 


momeee (Ff, F and the coefficient of eddy viscosity (Az ) 


y) 
Pemerated from (Vc). Auxiliary variables were defined at 
meoss {( X ) points. 

The model is not strictly designed to simulate any panee = 
tat ocean, but rather to be representative of the fundamental 
Physical characteristics of mass transport and western boundary 


Ere nuwOt an Ocean in the northern hemisphere as affected by 


Pic uESeOUS aAcCtLOM OF the ocean's large scale circulation. 
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Figure 1: The Staggered Grid with Corner Boundary. The solid line represents 
the physical boundary of the basin. 
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(@) Principal grid point (Ww, ¢,D, a(jpl) ) 
(o) Grid center (evan VG oc e eco Pio Fy ) 


(x) Grid average points (UAV, VAV, UAVE, VAVE, A(-{D])ove, eee yaa 
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The energy for this barotropic wind driven circulation 


model is represented by 


De 2ny 
a So = Jeane 5 ae Ty = 0, (III-2) 


pepoimmern Of westerly winds in the center half and easterly 


famasein the northern and southern quarters of the ocean 


f/(eceukigure 2). This leads to a stress curl term 
mo x 2nf Say: : 
pee a0 s a. A (aia 5) 


rime provides the actual forcing in the vorticity equation. 
faetsethe amplitude ©f the zonal component of the wind, and 
bemiice the north-south extent of the domain. 


lites verocity (Wav) written an numerical form is 


c 
i 


i,j GE eo ea) ; (esha SE] 
ee Cha : Vel (Yi-1,j ae) 


The boundary conditions ee (Una US ccm lene = EiNe Coastline 
were zero normal flow and zero slip. Zero normal flow was 
maeeamed by requiring the streamfunction (Y¥) on the left side 
Mmameaiation (ll-3) to be equal Zero on the ocean perimeter. 
immo lement the condition of zéro normal flow and zero slip 
mmmene terms On the right hand side of (11-3), the velocity 


Mmemactined aS zero on the coastline by defining the zonal 


boundaries 
ius Le SU OL, Gal = aa ay 
We ae Ol Vas =e. 35, (iets!) 


mer the meridional boundaries 
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Figure 2: Wind Stress Pattern 


oe 





wi, 3 Bae moe 33 ~ = SS55) 


a ee Mogi e = BS sae Cll to) 


and for the corners 
Cu,v)1 4 i Care 2 (u,v) 34 4 = (u,v)33 2 
Wey sAne Ge), 53 (VY) 34.34 7 (U,V)55.55 (111-7) 


Where the original 32 x 32 (u,v) grid field becomes a 


See oo) erid to include the boundary conditions. 


ieeNItE DIFFERENCE FORM OF THE VORTICITY EQUATION 
The finite difference form of (II1-3) was written 


y. y. ; 
oO g2y. . ( ete) 2 agaler)) 
ot eG we 2AX 


1 (V° Vv) G4) jlogs RO a 4 j 
yee 7 ) 


VON) ie ea + CVs VV) a 23 
oe 


ceealack ek! + (V-Vu) i#1,j+h) 


+ 
Qu) 


(eer We se eee Vlog 
— : ey 


F ° * ae F ° ° Ff © ° 4 FE ° ° 
omen re oe Yi+l,J) ae ae Y1 Jy] 
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aXe Gat fal Xe : Leer et eX ey 
: + (itl iis (1 
3 ae Sam aL yjJd (ues ) 
In (III-8) a is given by 


i 


2 = . . ° : : . = | es 
Sg a 1,3 -1,) * *4,5+1 * %4,5-1 ~ “Pajl- 
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The advection terms ( V-VU, etc.) were finite differenced 
feeme the Same metnod as “Haney [1974], and the friction 
Memimns were expressed differently for each of the two forms 
@epcddy ViISCOSIity aS discussed below. Readers are referred 
fmomene documented Computer program in Appendix A for further 
G@emallis, 

Minoo MmenCmilciI cel iemrs f1TrSt mecessary to determine 
three values of minimum A (AMIN), such that the maximum in 
Micmanalytic Stréamfunction (¥) 1s at x = d/2Z, ad, 2d, 
Mmoamceuwmyely. Since the grid cannot resolve the western 
Dbewliagary in the first case (where the western boundary width 
mee a large computational mode was expected with AMIN(1). 
Mmemsceond case was expected to be marginally stable with 
AMIN(2Z); and finally, with AMIN(3) it was expected that the 
teeeeri boundary would be clearly resolved, but the value of 
AMIN(3) would be so high as to unrealistically damp the 
mic~entor Océan solution. The respective difference equations 


for the three AMINS were: 


AMIN(1) = (B x “3 x d/2)3 
2n 
= v3 
AMIN(2) = (8 x 52 x d)3 
AMIN(3) = (g x A) x 2d)3 (Ciao) 
sal 


Mmese three values of AMIN were used for the three fundamental 
experiments where A was taken as constant, and also for the 
mominun A in the experiments with non linear coefficients of 


Pec e seCOcry sine (ii—-i5) and (11-24) with Ag = AMIN. 


Zo 





In determining the grid size and location for the non 
BiicHtiwcOctiMelcwpse or ncddy VauscoOsity, tne grid location of 
mie eecneTamne paraneters py: | or |D|) had to be taken into 
eomoideration. This led to a 34 x 34 grid for A(|Vc!) in 
Piemeasic O41 JCllstrTopny —GaSGade, and a 33 x 33 erid for A(| D)) 
Mmmeic Kiet G emersy cascade form. (See Figure 1). 

Pe Eiseh Op ivees Cadcm@Gasic 

LiMmnemoMcthOpiye cascade Case, im Order toO=generace 
wererrecLencts Of Non linear eddy viscosity, 1t is necessary 
Mmemoeoncravesa relatave vorticity fieldwand ito determine the 
Mer beciwtnO: VYOrtlicity at each time step. Several possible 
meomucauwes CX1St in developing a vorticity field, where in 


mm@meaises VOrLIicity would be defined on the 33 x 33 principal 


pao peints (°-). For this model € was defined in terms of 
Gite) by 
2 ea le oc 
OX oy 
a VAV i 41, 7 VAVS 3 ; UAVS 541 ; UAV; 
1,) d d (Clie 1e6)5) 


where UAV and VAV are defined as 


ee Jee! ee Nees Ve Ae 
ee i,j , Te EL) ; , 


Geet) 
cen) lec see 54g = «(1,433 


/memorder tO get an equivalent (u,v) on the principal grid 
fomees. Lt (UU wikareswritten in terms of ¥, then (III-10) 
meauUcesS tO 

Ci,j = wqe Viel, f+ ete i141 


Cute) 
Ws es 1 e +] —= Ys s 
se Ey, 





This method, along with three other schemes of Miyakoda [1962] 
HOrecalculating Vortieity, were utilized and compared. 
le 22a ite de eeOme heinuys (| Va mwase developed aising 


eateered Citierwemees on the 52 x 32) ¢rid centers ( °%); 
Sl GR rg = Stet este yg 
g 2 2 


a Ge ed oe es PP Oe 
2 2 


(ITI-13) 
WS Poe S85 IS By SE 
lye >) Seid co dititenence: orm 
AC|¥ol)4q 5 = AMIN x (1 + mx [Velj; 
IVo | max (III-14) 


iP = Poe |) = aero 
where Ajj is defined on the 32 x 32 grid centers. ines a pul 
of (ea Fae of 6 x 10'* was estimated from linear theory 
Oimak, 1950}, and m was considered an adjustable parameter. 
it it omdirmoneomnEorninet athe iri cLlon LoOnce sor 
ememrophy Cascade, cn a 32 x 32 grid, was written from (11-14) 
ame (11-15) as 


= L ;(Ai,j + Ait+1j (Ui+1,j - Yi,j) - 
. ° = — : > b 
FXG 5 q2 Lf a ; - ) ; J 15)) 
Ga Ld, Ly CUR gd Mata) 
o Ep Gini 2 Se ja ep 
: | 


= (A jt (Ut TL 5- 1) 
2 


4S: 





SM ae a | 
PYG 3 > aay Ed) (viel j - Yi) 
, ee 8) (7: , 5 = Wily 
Leas pre 
q2 EC A24 eee eh ye) 


- ol ieee | Vee) | (111-15) 
Mee. i SOY = tees SS 

Pow omelcdtmenac lige to) Tequares a S40 54 sited 
mae voivye|) in order to calculate the friction forces. Two 
meehods were utilized to generate a value of [Vc | and thereby 
mevalue Of eddy viscosity coefficients across the boundary. 
Mtestirst case was linear extrapolation in all four direéc- 
tions. For the western boundary the finite difference form 
of A(|vc|) was 


: SLOPE x d+ A, 5 (II1-16) 


a 
fe Se ae OO 
cpa the values of A immediately outside the other boundaries 
were determined in a similar manner. This method is physically 
Marymescitdtive because it gives a boundary of increasing ¢co- 
efficients in the direction of the western boundary. Along 
Omen other boundaries the gradient, and hence the boundary 
value of A(|VZ|) is flat. 
(ieee eoncmictmodmutntlized to, generate a value for 

gic) ) aeross the boundary was to extend the existing interior 


Meoundaiy Soutwands in all foursdirections, so that across the 


Wester Noundary, Lor example 
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in Nee 2s. OS Gey) 
This method is consistent with the antisymmetric velocity 
profile which accompanies the zero flow boundary conditions. 
In both cases the 34 x 34 corners were not used for IMGs [1a 
Hic lc eaemirncthon force. defined on a) 31 ox 15d 
Mie rior grid, was 


e . +r a e e . ° e « 
Epeeae Fi, j as (ees ; Dyagpel a Fyi-1,j ; fYa=1 3174 


See So ee 
tee Qa? ene 33 (III-18) 


Zoe tee tic Energy CaseadesCase 


Micon cciominedeteinamine the detormation fiedd 


Sueenme Llwidewas 1) toscalcullate the shearing deformation 


(Ds) 
eee AY Gs ee eee 
seg d GELS 1S) 
ime Sos 1S doc 5 BE 


fieere WAVY and VAV are defined in (III-il1), and 2) to calculate 


miemotrecching deformation (D,) 


De. . = UAVEi+1,; - UAVE],; _ VAVE: j+ VAVE3 , j 
1,j d sr 

(III -20) 

yee 1G OAS j=l. oS 

where UAVE and VAVE are defined as” 

UAVE ; } = Uj J Uj Lee jae a = eee 54 

J |) ae ee 

VAVE: =: = Vi j 7 Viel j ies egos 4,0: 
a 2 se 1... 34 (III -21) 


Then the square root of the squares of the values of 


ee 





SeosOrmation along and mormal to a streamline gave the value 


Serdeformation at Gach principal grid point (°): 


|D | eee) a Chey ey Gun 22) 


ea) 
i=1....33 j#41....33 


foots cdudatiom (Il1l-Z2Z5) im finite difference form was used to 


calculate values for A(| D]) Pesce deOilmnicnt te —eNeme yee ase adic 
A(|D|); a AMIN ome wets) (lili 23) 


| Pmax | 
ie nie Mee eS 


where AMIN and m are as described above, and the maximum 


deformation (D lmwdcecotimated. to be 6.22 x 10 ’.) Since 


max 


ai) Wasa 35 x S53 Eleld, the form of the friction force 
did not require additional boundary conditions for A(|D]). 
Miecwiumertcalwhorm ton the friction forces based on 


Peer) written from (I1-14) and (J1-15) was 


dle A a eee 


Slag ee [li Cae = Toa) 


2 

_) @Sieilaie es Set ily (CCL, 2 = Eat 
2 

1 A: - + A: : u: : u: : 

ed i Gsm! SUT) (Tea oe aes 
2 

= (Ai ,j a BSI 5 7} okay (ebb i Vs oe 

ae anes ee 


eee Aj 4- te gANS Wee wea \y = os 
Pyi jg © ge CCed ead) athe) © 45) 


- mei et : NLS sci gue ‘i wi-1,5)) 
2 
it ae ee re 
ike me oe) (toatl) 6 - Vi,j) 
Z. 
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ct. jee tel Yj > Yi, 5-1) 1 
<< a ; | (lel 2 a) 
The Smagorinsky form of the friction force for An |D| from 


(II-26) and (II-27) was 
Pxq jy (Aig * Ai,j-1)/2) [eq,5 * Pea, j-1)/2) 
Cee ete) See, 4 Pee? 24 

Cea eee cl} LCD. w UWsroi gd os 
25 AG 10/2) (sq-1 4-723} 

Fy; jzgfl i,j ee Se se Cn ere Dene 2. 
eC ete ee ate eye) ((Degey 5-4)7 21 
aD is at Aq-1,5)/21 pONea Y Deh sae ©. 


ea 2 Ada y1) 72) IDE. fol Deen gan’ ol: 


(eres) 
eee SiS 2a 3S 


ic wit Cud ii remenhGenOiinmomeurner clr or the friction 

moiecets the same as for the enstrophy cascade case (III-18). 
ee oOlLution Description 

Wsame sit) -S)ieas the basis for solution, a centered 
momemditterence scheme (Leapfrog) was used for all terms 
eweept the friction terms which were evaluated at the previous 
mimes step. With time steps of fourteen hours, equation (111-8) 
memeber raced tor 210 days, “lo Start the model and to prevent 
Zroruolmeseperation, the Eulér-Backward (Matsuno) time scheme 


Power ized tor tne first and vevery fifty time steps. In the 
29 





Kanctic energy cascade experiments the form of the friction 
force shown in the following results was in accordance with 
equation (I11-24). Experiments were also conducted using 
the Smagorinsky form (III-25) which gave similar results to 
(III-24) and therefore are not shown in this study. 
Pier cOUntToOn plase Ot the model jad to solye tne 

equation 

v2(Se) - F, = 0 Te 2G) 
for the tendency d¥/at. This was done using a direct 
meisson solver. This new technique was written by R. Sweet 
(1972], based on a method originated by Buneman [1969], and 
imeversed for this model by Professor F. Faulkner of the Naval 
Postgraduate school. The method is extremely accurate and 
made possible computer time savings of nearly an order of 


magnitude. 
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LV. RESULTS 


Pee Reo LAVEMEN!T OF PURPOSE 

ict PuENOse sO t Tre S chests Was to present. a Scheme 
Preeciy athe accliracy Oot the mumerical solution of a finite 
Prada ocean Circulation model would be improved by the intro- 
aaeronwo: Nom tinear lateral eddy viscosity coefficients: 
Ponioumin carinver Chapters, the gradient of relative 
WOmeteity Or the fluid deformation could be used as the 
maemueccilye Darameters to generate the coefficients of eddy 
fecosctmey. Ihe Solutions using Constant coefficients of 
lateral eddy viscosity will be compared with those using 
variable coefficients derived from enstrophy cascade (Av|Vz]) 


and kinetic energy cascade (Av |D|) respectively. 


Pee NALYTICAL CONSIDERATIONS 

(iicmiiiMitialwexpeCrImenis Investigated three cases of con- 
eeonit SoC im me fen (cmon wcdGay aVilscoslty =r irst, am accurate 
fieebytical solution for the streamfunction was made by means 
Sema separate model where grid spacing was 60 km. The ana- 
iteeed! Solution of the streamiunction in terms of the wind 
emeeess CUrl was developed by Munk [1950], who considered a 
ircar eddy Viscosity - 

Gry) = r(x) 8 = CURL (IV-1) 
where r = domain width, curl,t = the k component of the wind 


memess curl, and 


all 





Mix) = Ke~ ¥2kx ee, es + J] 


- cs Qap oe SOE = 2 1) 
in which X(x) = distance eastward from the western boundary, 
2 aes 
= 73 SE 
k = (6/A) ¥ and KK = ee > EF 


iMmemrcader 1S referred to Figure 3 which portrayss the 
Heeincical Streamtunetion ( pp) for the three cases where the 
Meena Of p occurs at d/Z, d and 2d, respectively. The 
Mmamec Values of A which permitted the above analytical situ- 
Berm nmCGCCurpmWwere tine three constant coerficients of eddy 


viscosity now examined in the numerical model. 


Sex PERIMENTS WITH CONSTANT COEFFICIENTS OF EDDY VISCOSITY 
lini ECON stalt COCta1Cient Of Eddy viscosity, 
Pee. t2 x 10° cm* sec *, which physically represents the 
mbeckior Ocean circulation most accurately, produces a larse 
Zeioude Computational oscillation which fills the entire 
Dee Oot the numerical model. Figure 4 shows the extent of 
miewoscillation produced in the t field by this relatively 
Moveiagnitude coefficient of eddy viscosity. Figure 7 shows 
a direct comparison of the analytical wy field and the numer- 
ical py field produced by A, at the latitude of maximum wind 
mmmess Curl. The reader 1S also referred to Table I which 
presents a tabular comparison of y field highlights as 
memerated by constant eddy viscosity coefficients. All of 
these features are in accordance with the study by Takano 


1975). 
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Neraee n= Ul oor x 10° cm* sec~*, an order of magnitude 
Harecr cociricient of eddy viscosity, was examined in the 
memerical model. Figure 5 Shows the circulation for constant 
Pweeeitgure 11, a graphical comparison of analytical and 
iimierxical Y field for Aj; and Table 1, a tabular comparison. 
fees, Used aS a constant, produced a marginally satisfactory 
molutton, where the western boundary current was rather well 
mepresemted, but a moderate computational oscillation was 
Still evident with the value of the maximum streamfunction 
Svomnmehner than the analytical solution. 

bymimercasing the coefficient yet another order of 
eweieeridemeorA. = (,5 X 10° cm* sec * (Figure 6, Figure 14, 
Baoelaple 1), the computational oscillation is negligible 
Wren only about 3% error in the numerical solution. The 
westerm boundary current placement was in accord with the 
Bimivymical case (Figure 3), but the solution in the ocean 


Miser Werewas unrealistically damped by the large viscosity. 


D. EXPERIMENTS WITH NON LINEAR COEFFICIENTS OF EDDY VISCOSITY 
Pee meotcCdtetidtsanerce #5 me saieble COCEFICIeEnt Of eddy 
Mmebee@ctey that can physically or numerically depict all the 
peyeectsS Of Fluid circulation both in the ocean interior and 
inthe western boundary. As can be observed from the results 
Pomerar, Che objective Of Wsing non linear eddy viscosity is 
mommave lOw COeLTICiIents Of Gddy viscosity in the ocean 
mieerlonr, and Increasing coefficients approaching the western 
moundary in order to resolve the western boundary current and 
Home vciiu tne Gevelopment of a computational mode in the 
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NOWmeuincar GCOCtiteremes that represent actual physical 
MmoGesscsein the ocean were developed in Chapters ii and 
muieeeein tS SHOWN be lOnmumac If Che limits of the range of 
these coefficients are properly chosen the desired objective 
Mo Acire ved . 

Mme aT Ss anioMmmieiedwiecxperimentsS, ranges of coeffi - 
cients of eddy viscosity were chosen with the minimum 


coefficient equal to AMIN(1) = A,, and the range of variation 


] 
Ou tie CoOecLfiCients was governed by the adjustable parameter 
iiimmocin, tl-lhole for enstrophy cascade and (11-25) for 
[mimecuLGmenerey Gascade. Of course, in this case where AMIN 
is the smallest, the greatest range of m was required to 
Monet iy resolve the western boundary and to prevent an un- 
acceptable computational oscillation from developing in the 
Pmeron es EXPeriments were conductéd with m varying from 
Mimeome00 with results that are noted below. It should be 
ioeedmenrat the value for m iS not précisely a direct multi- 
Pelctmiror the range of A due to the non linear effect of 


ie cli loi) and iiss) The actual 


lVomax| ren 


fees OL A/AMIN were printed out in the experiments and the 
range of A/AMIN appears in Table II for the most significant 
eespecriments. 

In the enstrophy cascade (Av |Vc]) experiment for 
A, = 0.12 x 10° and m = 10 and using the symmetrical boundary 
moMmcalilons or A as imdicated wm (ilil-17), the actual range 
esi, was 1.0 to 7.0. The resulting 7 field can be seen 


Mime rsn. Usame the lamegr extrapolation shown in 
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Symmetric Boundary Conditions for Av|Vvz| 


NMUMel Ve Mies oemitlOlowouw)s for A, 


repeure Gil: 
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weit 6) tO Obtain Avacross the boundary, the solution shown 
in Figure 8B was obtained. In this case the range of A/A, 
ee tO /.9. Slt etSem@ted ere that an the extrapolation 
boundary condition the maximum A occurs outside the western 
boundary and therefore is artificially derived; and that the 
A/A, minimum values indicated in these experiments are not 
meeuabiy the lowest ratio obtained, but a value representative 
Geeche A obtained in the interior ocean. In the kinetic energy 
wascadesexperiment (Aw |pi) for this case, the range of A/A, 
eee co 1058) with results very similar to the enstrophy 
@iscade Experiments (Figure 8C). As can be readily Seen) oud 
wtiemireures for the above three experiments with AMIN = A, 
aman - 10, a heavy computational oscillation is still very 
Imei cxistence, although a definite improvement over 
Pemaeconotant (Figure 4) is apparent. 

memcaneaiso be noted here and in the following experi- 
ments, that although the method used in deriving the non 
imledmecoctficients of eddy viscosity, namely A» |Vc| and 
Au|D|, show somewhat different characteristics in the ocean 
MeeeiideLon pattern, the results were analogous enough that 
Baemiiain purpose of this thesis could have been accomplished 
taeieettner method and either boundary conditions for A in 
the case of Av|Vc|. In addition, several tests using various 
Sompinatwons of the Laplacian given by (111-12) and the usual 
S5-point Laplacian were made. It was found that the method of 
feminine the relative vorticity field was of no consequence 


Mimereshesmics gcibeved, and therefore this paper does not 
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Teer ant any curener discussion Of this aspect. In all the 
above non linear experiments with A,, in addition to the 
imeonines listed and Table Ils the reader is referred to the 
peaphical representation of the steady state maximum analyt- 
ical and numerical streamfunction field in Figure 7. 

In the next set of experiments with AMIN(1), m is in- 
m~eascd to 100 civing a two order of magnitude range for the 
Seerricicnts of eddy yiscosity. The results are shown in 
lmpoeeos IA, 9B, 9C, 7, and Table JI. In the cases where 
Av|V¢| with the symmetric boundary conditions for A, A/A, 
mamecd trom 1.2 to 41.7, and for the linear extrapolation 
boundary conditions A/A, ranged from 1.3 to 62.2. For Av|D| 
miemvalues of the coefficients ranged from 1.5 to 70.1. In 
Meese orOuUp Of experiments the solution of streamfunction 
moroved preatly, with the numerical Da Vet beTe cede Wie Okoae 


Somme Of the analytical wp in all cases. However, the higher 


max 
we OL A at the western boundary resulted in a tendency for 


Y 


mee co move eastward, especially in the case of Av|D|. 
mieneasing the range of A) another order of magnitude, 
Minmiext CxXperiments examined m = J000. These results are 
Mamet Pigures 10A, LOB, 10C, 7, and Table JI. In all cases 
the computational oscillation was negligible, the interior 
solution was not strongly damped, and the western boundary 
fmeappen was well defined. The variation of A/A, was 1.5 to 

ies for Av| Ve | (Symmetricehoundary  Genaditions), 1.6 to 23/7 


for Av|Vc| (extrapolated boundary condition), and 5.0 to 313 


for Av[D|. Again, however, there was a tendency for p,,, to 
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Brounen OG: wWumerical Solutions, ore) for A, and 
m = 1000, Ar|D| 
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Nene castward to 3c in the case of Ay|Di. Due to this, @and 
Gicesligchtiy higher minimum Coefiiciéents in the ocean interior 
momecne An |D| solution, the Au|Vcg| solution appears to be 
memerally preterred to Avi Di) This is probably because the 
method based on two dimensional turbulence is somwhat more 
womcative to the space scalefof motion. 

Bn the next sets of experiments, AMIN was increased to 
mee- 0.95 x 10° cm* sec ? and m was examaned for values ranging 
maom 20 to 100. These results are shown in Figures#ll, 1Z2A, 
meee tS, and also Table III. With a higher initial AMIN and 
fmecwitineg higher viscous solution in the ocean interior, the 
computational oscillation was Suppressed and a completely 
Saerstactory solution was attained by m = 20. The range of 
meoues for A/A, and y,,, are given in Table 111, and Figure 11 
meeees a Graphical representation of results for the non linear 
eeaecriments with A,. Inereasing m to 100 had the undesirable 
effect of moving the western boundary eastward to 3d in the 
Base Ot Aw|D| . 

Non linear experiments with AMIN = A, = 7.5 x 10° cm’ sec } 
produced no enhancing results. The solution was exceptionally 
viscous, and the interior ocean was already overdamped with 
Poi? moh minimum coetficient of Cday viscosity. Increasing 


m to 10 resulted in moving the wy Home| | Castward £06 


max 
Sd . Rememeto Taplem ll toregeprescwumenve values of A/ A. 
mae ,, and to Figure 14 for the streamfunction field for 


A~w|D| with m = 10. 
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Pieeauchnor 15 O02 the Opinion tides tie wecpcrincitiew tem 
AMIN = A, (or perhaps a slightly higher minimum coefficient) 
and m ranging somewhere between 100 and 1000 (depending on 
minimum acceptable computational oscillation) gives the best 
mild Of non linear lateral eddy viscosaty coefficients for 
iomemoe Utilization in more sophisticated finite difference 
Pecomecirculation models; [t appears from the experiments 
that the only drawback is that these higher values of m 
Pmeauce a western boundary width which is nearly 2d. 

An accurate eee and mumerireal cdepretromvoe bots tiie 
Ocean interior and the western boundary with no computational 
eoenllation was achieved by using either of the two forms of 
iemelinear eddy viscosity. These were achieved with the 
heiimum Coefficients approximately equal to A, in the interior 
Seoam, increasing approximately two orders of magnitude in 


the more dynamic flow in the western boundary region. 
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V. CONCLUSIONS 


This numerical, mode hasebeen slecess htm best Tne «a 
Fenemeawienreby thewacctiracy o: the nNimmeriecal solutton of a 
mimite grid ocean circulation model can be improved by the 
introduction of non linear lateral eddy viscosity coefficients. 
Non linear coefficients, properly generated to represent 
Meulcdmpny sical processes that may be occurring in the ocean, 
alow the use of low Cocificilents in the Ineerior ocean solu- 
mioneamd high coefficients in the higher circulation density 
mene ewestern boundary Current. This distribution of eddy 
femecOstty 1S SUffICient to prevent the formation of a compu- 
fGational oscillation which would occur if the low value were 
used throughout vhs domain. 

Two methods, namely enstrophy cascade (Av |Vz|) and kinetic 
energy cascade (Av |D]), were investigated and presented which 
teeallow the researcher to generate non linear coefficients 
in other models. From experimentation with this barotropic 
ieael, minimum coefficients and corresponding ranges of 
G@ecetficient magnitudes are recommended for initial investiga- 
tion in more sophisticated baroclinic ocean and coupled 


weean-atmosphere finite difference models. 
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